Building Systems with Roto-Translations

For building systems programatically, PyBigDFT offers some helper routines based on the concept of rototranslations. Using these features, we can easily align molecules in space.

Rotations and Translations

Let’s begin by showing the basics of rotating and translating molecules.

from BigDFT.IO import XYZReader
from BigDFT.Systems import System
from BigDFT.Fragments import Fragment

sys = System()
sys["WAT:0"] = Fragment()
with XYZReader("H2O") as ifile:
    for at in ifile:
        sys["WAT:0"] += Fragment([at])

Translation. Note that the units are atomic units.

from copy import deepcopy
sys["WAT:1"] = deepcopy(sys["WAT:0"])
sys["WAT:1"].translate([-10, 0, 0])

from BigDFT.Visualization import InlineVisualizer
viz = InlineVisualizer(400,300)
viz.display_system(sys)

You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol

Rotation. We can pick our units of degrees or angstroems.

sys["WAT:1"].rotate(x=90, units="degrees")
viz = InlineVisualizer(400,300)
viz.display_system(sys)

You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol

Fragment Interpolation

We can interpolate between two fragments if we wish to examine some intermediary states.

from BigDFT.Fragments import interpolate_fragments
steps = interpolate_fragments(sys["WAT:0"], sys["WAT:1"], steps=3)

sys2 = System()
for i, s in enumerate(steps):
    sys2["WAT:"+str(i)] = s

viz = InlineVisualizer(400,300)
viz.display_system(sys2)

You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol

Lining Up Fragments

Sometimes we have two fragments which are made of the same atoms, but one geometry is perturbed from another. We might wish to try lining up those fragments so we can figure out how big that pertubation is.

from BigDFT.Fragments import RotoTranslation
rtsys = System()
rtsys["WAT:0"] = deepcopy(sys["WAT:0"])
rtsys["WAT:1"] = deepcopy(sys["WAT:0"])
rtsys["WAT:1"].translate([-5.0, 0, 0])
pos = rtsys["WAT:0"][0].get_position()
pos[1] -= 0.75
rtsys["WAT:0"][0].set_position(pos)

viz = InlineVisualizer(400,300)
viz.display_system(rtsys)

You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol

Using the lineup freature, we can find the best matching between these two systems.

from BigDFT.Fragments import lineup_fragment
rtsys["WAT:0"] = lineup_fragment(rtsys["WAT:0"])
rtsys["WAT:1"] = lineup_fragment(rtsys["WAT:1"])

viz = InlineVisualizer(400,300)
viz.display_system(rtsys)

You appear to be running in JupyterLab (or JavaScript failed to load for some other reason). You need to install the 3dmol extension:
jupyter labextension install jupyterlab_3dmol

Now compute the RMSD.

from numpy import array
from numpy.linalg import norm

rmsd = 0
for at1, at2 in zip(rtsys["WAT:0"], rtsys["WAT:1"]):
    rmsd += norm(array(at1.get_position("angstroem")) - array(at2.get_position("angstroem")))

print(rmsd)

3.681664737333418